Qual é a alturas média das årvores em uma floresta?
Departamento de Ecologia IB-USP
Qual é a alturas média das årvores em uma floresta?
Fazer afirmaçÔes sobre as caracterĂsticas de uma população, com base em informaçÔes dadas por amostras.
População Ă© qualquer conjunto de elementos que vocĂȘ estĂĄ observando
Amostra é qualquer subconjunto desta população.
## [1] 22.72789
## [1] 32.3065
## [1] 5.683881
A variĂĄvel aleatĂłria X Ă© uma variĂĄvel que tem um valor Ășnico (determinado aleatoriamente) para cada resultado de um experimento (lembre-se esse termo veio da teoria de probabilidade).
Exemplos:
Altura das ĂĄrvores de uma floresta (!!)
NĂșmero de presas capturadas por um predador em um determinado dia;
Comprimento de um peixe adulto selecionado aleatoriamente.
As variĂĄveis aleatĂłrias podem ser discretas ou contĂnuas.
Função que define as probabilidades P(X) para cada valor possĂvel da variĂĄvel X.
Distribuição de Probabilidades.
Discreta: nĂșmero de sucessos em n tentativas . p - prob sucesso.
NĂșmero sementes predadas em cada experimento contendo 20 semente.
NĂșmero de animais infectados em relação ao nĂșmero de capturados em cada local
Discreta: Probabilidade de uma sĂ©rie de eventos ocorrem em um perĂodo fixo de tempo, ĂĄrea, volume, quadrante, etc.
Dados de contagem
Abundùncia de espécies em um local
NĂșmero de indivĂduos capturados por tempo
dpois(5,10)â Metade (50%) dos dados estĂĄ acima (e abaixo) da mĂ©dia
â 68% dados estĂĄ dentro de 1 desvio padrĂŁo da mĂ©dia
â 95% estĂĄ dentro de 2 desvios padrĂ”es da mĂ©dia
â Virtualmente todos os valores estĂŁo dentro de 3 desvios padrĂ”es da mĂ©dia
## [1] 0.005856729
"Toda soma de variĂĄveis aleatĂłrias independentes de mĂ©dia finita e variĂąncia limitada Ă© aproximadamente Normal, desde que o nĂșmero de termos da soma seja suficientemente grande."
Independentemente do tipo de distribuição da população, na medida em que o tamanho da amostra aumenta, a distribuição das médias amostrais tende a uma distribuição Normal.
levels(arv$sp)
## [1] "sp1" "sp2"
Serå que a altura média de uma espécie nessa floresta é diferente da altura média da outra espécie?
Ou seja, serå que as espécies possuem uma mesma distribuição de probabilidades e as diferenças encontradas são devido ao acaso ou serå que cada espécie é uma variåvel aleatória com médias diferentes?
mean(arv$alt[arv$sp=="sp2"]) - mean(arv$alt[arv$sp=="sp1"])
## [1] 5.109466
Posso afirmar que as médias são diferentes?
H0 - não hå diferença entre as alturas das duas espécies. A média da altura da espécie 1 é igual à média da altura da espécie 2 e a diferença observada se deve ao acaso.
H1 - hå diferença entre as alturas das duas espécies. A média da altura da espécie 1 é diferente da média da altura da espécie 2.
AFIRMAĂĂO: nĂŁo hĂĄ diferença entre as alturas das espĂ©cies,
VERDADE: hå diferenças.
AFIRMAĂĂO: hĂĄ diferença entre as alturas das espĂ©cies,
VERDADE: não hå diferenças.
ERRO DO TIPO I (α): rejeitar a hipótese nula dada que ela é verdadeira.
ERRO DO TIPO II (B): aceitar a hipĂłtese nula dado que ela Ă© falsa.
Contrastar a probabilidade de que a diferença que vocĂȘ encontrou nas mĂ©dias das alturas das duas espĂ©cies pode se considerada devido ao acaso ou nĂŁo.
Qual é a probabilidade de que as duas amostras pertencem à mesma população de medidas?
Podemos ou não rejeitar a hipótese nula (de que a diferença é ao acaso), baseado em alguma probabilidade de estarmos cometendo algum erro (α).
"O nĂșmero de observaçÔes independentes menos o nĂșmero de parĂąmetros estimados."
Quantas observaçÔes independentes foram usadas para calcular a média e a variùncia dos dados sob a hipótese nula?
100 ĂĄrvores
2 parĂąmetros
98 graus de liberdade
NĂŁo entraremos em detalhes de como se calcula a estatĂstica t e como Ă© a sua distribuição de probabilidades (recomendo ver na wikipedia). O importante a saber aqui Ă© que tendo as duas amostras modeladas como uma distribuição normal, calculamos a estatĂstica t e com um certo valor de graus de liberdade, nĂłs podemos recorrer Ă s tabelas da distribuição t para ver qual a probabilidade de se obter aquele valor que obtivemos dada que a hipĂłtese nula Ă© verdadeira. Se esta probabilidade for muito pequena, a chance de estarmos caindo no erro do tipo 1 Ă© pequena.
No R, a função usada para fazer um teste t é t.test.
t.test(arv$alt~arv$sp, var.equal=T)
## ## Two Sample t-test ## ## data: arv$alt by arv$sp ## t = -5.0125, df = 98, p-value = 2.387e-06 ## alternative hypothesis: true difference in means is not equal to 0 ## 95 percent confidence interval: ## -7.132311 -3.086621 ## sample estimates: ## mean in group sp1 mean in group sp2 ## 20.17316 25.28262
#o argumento var.equa=T diz que as variùncias das duas espécies são as mesmas # podemos fazer o teste com a premissa de que as variùncias não sejam iguais, # usando var.equal=F, isso vai "consumir" graus de liberdade.
Vamos interpretar este resultado: temos o valor de t calculado e os graus de liberade (jå calculado pela gente antes). Com estes 2 valores, o que o R faz é calcular a probabilidade de que a hipótese nula seja verdadeira. Esse é o valor-p, que neste exemplo deu um valor beem pequeno. Vemos então que podemos rejeitar H0, com uma probabilidade muito baixa de estarmos incorrendo no erro de tipo 1, e desta forma acreditamos que as espécies de årvores tenham distribuiçÔes de alturas diferentes (médias diferentes).
E se tivĂ©ssemos na verdade 3 espĂ©cies de ĂĄrvores para compararmos as alturas? Acabamos de ver que o teste t serve apenas para comparar duas amostras. O que chamamos de modelos lineares Ă© uma grande famĂlia de modelos estatĂsticos que modelam uma relação linear entre a variĂĄvel dependente (Y) e a(s) variĂĄvel(is) independente(s). Nesta famĂlia encontra-se o teste t que acabamos de ver, a AnĂĄlise de VariĂąnci (ANOVA), a regressĂŁo, etc.
Chamamos de Modelos Lineare Gerais (LMs), aqueles que possuem como premissa a variĂĄvel dependente Y como de distribuição Normal, e de Modelos Lineares Genearlizados (GLMs) os que assumem outros tipos de distribuição para a variĂĄvel dependente (como a Binomial e Poisson). No R, a função que usamos para estes modelos sĂŁo lm e glm (hĂĄ alguns mais especĂficos como o prĂłrpio t.test e o aov para ANOVA). O R lida com todos os modelos lineares de maneira similar, e um detalhe importante Ă© que todas as anĂĄlises desta famĂlia devem ser salvas em um objeto, para seu resultado ser apresentado apĂłs o comando summary.
Ă importante lembrar que neste roteiro estamos lidando com modelos paramĂ©tricos, ou seja, que modela a variĂĄvel de interesse como uma variĂĄvel aleatĂłria pertencente a uma certa distribuição de probabilidades. Existem outras abordagens estatĂsticas nĂŁo-paramĂ©tricas que fazem testes que nĂŁo tem como premissa a distribuição de probabilidades. Veja um exemplo de teste t nĂŁo paramĂ©trico neste roteiro.
Mesmo fazendo parte de uma mesma famĂlia de modelos, nas sessĂ”es seguintes vamos separar os modelos lineares pelos nomes tradicionais das anĂĄlises: anĂĄlise de variĂąncia, regressĂŁo linear simples e anĂĄlise de covariĂąncia. Depois, veremos alguns exemplos de GLMs aplicĂĄveis a dados ecolĂłgicos.
Vamos entĂŁo incluir no nosso exemplo inicial uma terceira espĂ©cie amostrada (e vamos diminuir o nĂșmero de amostras por espĂ©cie para ficar mais fĂĄcil de fazer os cĂĄlculos passo a passo):
alt.sp3 <- c(36.90076761,33.23131727,32.82767311,24.93410577,31.87626329,
28.76824248,26.6436144,31.41238398,27.65383929,31.70425719,
33.60752445,28.30263354,29.32210674,23.4790054,22.20793046,
18.28797293,23.39044736,28.75820917,29.94344703,34.47430235,
32.09820862,26.99845592,27.64858951,31.90557306,33.63236368,
25.35966179,35.00885172,32.89125598,18.98802229,30.8516906,
22.04191398,35.37232889,35.67150408,34.30257037,24.43038475,
20.6274663,31.15717916,23.27023231,36.20568545,34.33068448,
28.40769054,31.00578874,29.50989158,31.69915991,37.25975797,
30.1279997,38.623371,35.36472829,29.88304259,31.00145491)
sp3 <- data.frame(alt=alt.sp3, sp=rep("sp3",50))
arv2 <- rbind(arv,sp3)
amostras <- seq(1,150,by=5)
arv3 <- arv2[ amostras, ]
Neste exemplo lidaremos com um modelo simples de anova chamado one way ANOVA, porque lida apenas com uma variĂĄvel independente categĂłrica (fator).
A ANOVA testarå a hipótese nula de que as médias das alturas das 3 espécies não diferem.
Faremos agora a ANOVA "na unha", para entendermos melhor cada passo da anålise e a interpretação de seus resultados para depois usarmos a função do R que faz esta anålise diretamente.
Vamos calcular os valores da tabela de ANOVA. Começando com os desvios quadrĂĄticos, ou seja, quanto os dados desviam da mĂ©dia (idĂ©ia parecida com a variĂąncia). O ponto importante Ă© que essa variação Ă© aditiva e portanto, pode se decomposta. A variação total Ă© decomposta em variação relacionada ao tratamento (entre grupos), no nosso caso Ă s espĂ©cies, e uma variação interna dos grupos (chamada de erro). A estatĂstica F Ă© a razĂŁo entre essas variaĂ”es apĂłs dividir cada uma delas pelos seus respectivos graus de liberdade. Complicou? Vamos fazer os cĂĄlculos e ver os grĂĄficos para ver se entendemos melhor:
A tabela de ANOVA que vamos construir Ă© essa:
Primeiro vamos mudar a forma do nosso data frame para facilitar os comandos
arv4 <- data.frame(sp1=arv3$alt[arv3$sp=="sp1"],
sp2=arv3$alt[arv3$sp=="sp2"],
sp3=arv3$alt[arv3$sp=="sp3"])
arv4
## sp1 sp2 sp3 ## 1 21.282000 19.14508 36.90077 ## 2 19.875790 18.66252 28.76824 ## 3 23.562570 25.10581 33.60752 ## 4 19.313744 18.10266 18.28797 ## 5 26.308616 21.48084 32.09821 ## 6 8.301795 28.08559 25.35966 ## 7 27.192502 20.47805 22.04191 ## 8 15.941866 29.22927 20.62747 ## 9 20.596624 19.32911 28.40769 ## 10 22.048945 18.21446 30.12800
boxplot(arv4)
Vamos calcular a média geral, as médias para cada espécie, e as diferenças entre a média geral e a média para cada espécie, para chegar ao valor da soma dos quadrados total:
media.geral <- mean(arv3$alt) media.geral
## [1] 23.28284
medias.sps <- apply(arv4, 2, mean) medias.sps
## sp1 sp2 sp3 ## 20.44245 21.78334 27.62274
dif.geral <- arv4 - media.geral ss.especies <- dif.geral^2 ss.especies
## sp1 sp2 sp3 ## [1,] 4.00337202 17.121094 185.447876 ## [2,] 11.60801175 21.347420 30.089610 ## [3,] 0.07824755 3.323227 106.599051 ## [4,] 15.75374565 26.834280 24.948725 ## [5,] 9.15530094 3.247214 77.710675 ## [6,] 224.43179074 23.066347 4.313177 ## [7,] 15.28543337 7.866884 1.539904 ## [8,] 53.88994194 35.359990 7.051024 ## [9,] 7.21577099 15.631997 26.264064 ## [10,] 1.52250306 25.688498 46.856173
ss.total <- sum(ss.especies) ss.total
## [1] 1033.251
Para calcular os desvios quadrĂĄticos totais, nĂłs subtraĂmos cada altura das ĂĄrvores pela mĂ©dia geral e elevamos ao quadrado. Veja o grĂĄfico:
vetor.cor <- rep(1:3, each=10) #vetor de cores
plot(x = c(1:30), y = arv3$alt, ylim=c(10,40),
pch=(rep(c(15,16,17), each=10)),
col=vetor.cor,
ylab="Variåvel Resposta", xlab="ObservaçÔes",
main="Variação total")
for(i in 1:30)
{
lines(c(i,i),c(arv3$alt[i],mean(arv3$alt)),col=vetor.cor[i])
}
abline(h=media.geral)
Agora vamos fazer a somatĂłria dos desvios quadrĂĄticos dentro de cada grupo (ss.intra).
O grĂĄfico para entender esse cĂĄlculo:
vetor.medias<-rep(medias.sps, each=10)
plot(c(1:30), arv3$alt, ylim=c(10,40),
pch=(rep(c(15,16,17),each=10)),
col=vetor.cor,
main="Variação Intra Grupos",
ylab="Variåvel Resposta", xlab="ObservaçÔes")
for(i in 1:30)
{
lines(c(i,i),c(vetor.medias[i],arv3$alt[i]),col=vetor.cor[i])
}
lines(c(1,10),c(medias.sps[1],medias.sps[1]),col=1)
lines(c(11,20),c(medias.sps[2],medias.sps[2]),col=2)
lines(c(21,30),c(medias.sps[3],medias.sps[3]),col=3)
CĂĄlculo dos valores:
medias.sps
## sp1 sp2 sp3 ## 20.44245 21.78334 27.62274
ss.sp1=sum((arv4$sp1-medias.sps["sp1"])^2) ss.sp1
## [1] 262.2655
ss.sp2=sum((arv4$sp2-medias.sps["sp2"])^2) ss.sp2
## [1] 157.0018
ss.sp3=sum((arv4$sp3-medias.sps["sp3"])^2) ss.sp3
## [1] 322.4728
ss.intra=ss.sp1+ss.sp2+ss.sp3 ss.intra
## [1] 741.7401
A soma dos quadrados entre grupos:
plot(c(1:30), vetor.medias, ylim=c(10,40),
pch=(rep(c(15,16,17),each=10)),
col=vetor.cor,
main="Variação Entre Grupos",
ylab="Variåvel Resposta", xlab="ObservaçÔes")
for(i in 1:30)
{
lines(c(i,i),c(vetor.medias[i],mean(vetor.medias)),col=vetor.cor[i])
}
abline(h=media.geral)
points(c(1:30),arv3$alt, ylim=c(10,50),
pch=(rep(c(0,1,2),each=10)), col=vetor.cor, cex=0.5)
medias.sps
## sp1 sp2 sp3 ## 20.44245 21.78334 27.62274
media.geral
## [1] 23.28284
ss.entre=10*sum((medias.sps-media.geral)^2) ss.entre
## [1] 291.5112
#conferindo os cĂĄlculos ss.intra+ss.entre
## [1] 1033.251
ss.total
## [1] 1033.251
CĂĄlculo do F
# Desvios médios ms.entre=ss.entre/2 ms.intra=ss.intra/27 ms.entre
## [1] 145.7556
ms.intra
## [1] 27.47186
#F - razĂŁo das variĂąncias F.sps=ms.entre/ms.intra F.sps
## [1] 5.305634
Vamos ver na distribuição F qual a probabilidade de termos econtrado o valor F.sps ao acaso:
curve(expr=df(x, 2,27),
main="Distribuição F de Fisher (df=2,27)",
xlab="Valor F",
ylab="Densidade ProbabilĂstica (df)",
xlim=c(2,12))
abline(v=F.sps, col="red")
abline(h=0, lty=2)
xf=seq(F.sps, 12, 0.01)
ydf=df(xf, 2, 27)
polygon(c(F.sps,xf),c(0,ydf),col="red")
text(x= 7,y=0.08,paste("pf(x) =",
round(pf(F.sps,2,27,lower.tail=F),4)),
cex=0.8, col="red")
CĂĄlculo do P
p.sps=pf(F.sps, 2, 27, lower.tail=FALSE) p.sps
## [1] 0.01139248
Para mais informaçÔes sobre o F e as comparaçÔes com o teste t, veja esse roteiro.
Colocando os dados calculados em nossa tabela de anova
A função que constrĂłi esta tabela de ANOVA e faz o teste estatĂstico Ă© a aov. Vamos comparar:
anova.sps <- aov(alt~sp, data=arv3) summary(anova.sps)
## Df Sum Sq Mean Sq F value Pr(>F) ## sp 2 291.5 145.76 5.306 0.0114 * ## Residuals 27 741.7 27.47 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Agora nĂłs sabemos de ondem vem os nĂșmeros nesta tabela e como foi calculado o P, certo?
Onde estão as diferenças?
veja os grĂĄficos
faça teste à posteriori - Tukey HDS
TukeyHSD(anova.sps)
## Tukey multiple comparisons of means ## 95% family-wise confidence level ## ## Fit: aov(formula = alt ~ sp, data = arv3) ## ## $sp ## diff lwr upr p adj ## sp2-sp1 1.340893 -4.4708805 7.152667 0.8360276 ## sp3-sp1 7.180300 1.3685259 12.992073 0.0132022 ## sp3-sp2 5.839406 0.0276327 11.651180 0.0487468
#revendo grĂĄfico dos dados
stripchart( arv3$alt~arv3$sp,
vertical = TRUE, pch = 16,
col = c("black", "red", "green"))
#boxplot boxplot(arv4)
lm.anova.sp <- lm(alt~sp, data=arv3) summary.aov(lm.anova.sp) #mesmo que o summary do aov
## Df Sum Sq Mean Sq F value Pr(>F) ## sp 2 291.5 145.76 5.306 0.0114 * ## Residuals 27 741.7 27.47 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# mostra o resultado de forma um pouco diferente, consegue entender? summary(lm.anova.sp)
## ## Call: ## lm(formula = alt ~ sp, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -12.1407 -3.0002 -0.0742 3.2719 9.2780 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 20.442 1.657 12.334 1.32e-12 *** ## spsp2 1.341 2.344 0.572 0.57202 ## spsp3 7.180 2.344 3.063 0.00492 ** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 5.241 on 27 degrees of freedom ## Multiple R-squared: 0.2821, Adjusted R-squared: 0.229 ## F-statistic: 5.306 on 2 and 27 DF, p-value: 0.01139
plot(arv3$alt ~ nutri)
alt.nutr <- lm(alt ~ nutri, data = arv3) summary(alt.nutr)
## ## Call: ## lm(formula = alt ~ nutri, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -6.7176 -1.7171 0.0088 1.4124 9.2072 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 14.38711 1.05832 13.594 7.44e-14 *** ## nutri 0.58006 0.05971 9.714 1.82e-10 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.906 on 28 degrees of freedom ## Multiple R-squared: 0.7712, Adjusted R-squared: 0.763 ## F-statistic: 94.37 on 1 and 28 DF, p-value: 1.818e-10
plot(arv3$alt ~ nutri) abline(alt.nutr, col = "red") # reta da regressĂŁo #reta da hipĂłtese nula de ausĂȘncia de efeito dos nutrientes abline(h = mean(arv3$alt), col = "blue", lty = 2)
Y continuo
X1 contĂnuo
X2 categĂłrico
plot(alt~nutri, data=arv3, col=arv3$sp)
alt.sp.nutr <- lm(alt ~ nutri + sp, data = arv3) summary(alt.sp.nutr)
## ## Call: ## lm(formula = alt ~ nutri + sp, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -5.6682 -0.9554 -0.1585 1.3581 7.6134 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 13.3971 1.1034 12.142 3.23e-12 *** ## nutri 0.5256 0.0561 9.369 8.09e-10 *** ## spsp2 1.6421 1.1423 1.437 0.16251 ## spsp3 3.8326 1.1965 3.203 0.00357 ** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.553 on 26 degrees of freedom ## Multiple R-squared: 0.836, Adjusted R-squared: 0.817 ## F-statistic: 44.16 on 3 and 26 DF, p-value: 2.401e-10
alt.sp.nutr2 <- lm(alt ~ nutri + sp + nutri:sp, data=arv3) alt.sp.nutr3 <- lm(alt ~ nutri*sp, data=arv3) summary(alt.sp.nutr3)
## ## Call: ## lm(formula = alt ~ nutri * sp, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -5.6724 -0.9444 -0.3313 1.0416 7.0702 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 13.40167 1.42727 9.390 1.66e-09 *** ## nutri 0.52527 0.08906 5.898 4.38e-06 *** ## spsp2 2.93717 1.97080 1.490 0.149 ## spsp3 0.43636 2.75666 0.158 0.876 ## nutri:spsp2 -0.10095 0.12423 -0.813 0.424 ## nutri:spsp3 0.17187 0.14350 1.198 0.243 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.474 on 24 degrees of freedom ## Multiple R-squared: 0.8578, Adjusted R-squared: 0.8282 ## F-statistic: 28.96 on 5 and 24 DF, p-value: 2.003e-09
summary.aov(alt.sp.nutr3)
## Df Sum Sq Mean Sq F value Pr(>F) ## nutri 1 796.8 796.8 130.179 3.53e-11 *** ## sp 2 66.9 33.5 5.467 0.0111 * ## nutri:sp 2 22.6 11.3 1.846 0.1796 ## Residuals 24 146.9 6.1 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Outras distribuiçÔes de probabilidades além da Normal
VariĂąncias nĂŁo homogĂȘneas
Mais comuns:
Poisson â Ășteis para dados de contagem
Binomial â Ășteis para dados com proporçÔes, ou presença/ausĂȘncia
Uma hipótese sobre a distribuição da variåvel resposta Yi.
Especificação do modelo (variåveis X).
3.Função de ligação entre a média Yi e as variåveis X. Função que lineariza a regressão
| Distribuição | link function |
|---|---|
| normal | identity |
| poisson | log |
| binomial | logit |
| Gamma | reciprocal |
Y - NĂșmero de atropelamentos de anuros em uma estrada
X - distĂąncia a um Parque Nacional.
plot( atrop$TOT.N ~ atrop$D.PARK,
xlab = " DistĂąncia ao parque",
ylab = "NĂșmero de atropelamentos")
mod <- glm(TOT.N ~ D.PARK, data = atrop,
family = poisson)
summary(mod)
## ## Call: ## glm(formula = TOT.N ~ D.PARK, family = poisson, data = atrop) ## ## Deviance Residuals: ## Min 1Q Median 3Q Max ## -8.1100 -1.6950 -0.4708 1.4206 7.3337 ## ## Coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) 4.316e+00 4.322e-02 99.87 <2e-16 *** ## D.PARK -1.059e-04 4.387e-06 -24.13 <2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## (Dispersion parameter for poisson family taken to be 1) ## ## Null deviance: 1071.4 on 51 degrees of freedom ## Residual deviance: 390.9 on 50 degrees of freedom ## AIC: 634.29 ## ## Number of Fisher Scoring iterations: 4
anova(mod, test = "Chisq")
## Analysis of Deviance Table ## ## Model: poisson, link: log ## ## Response: TOT.N ## ## Terms added sequentially (first to last) ## ## ## Df Deviance Resid. Df Resid. Dev Pr(>Chi) ## NULL 51 1071.4 ## D.PARK 1 680.55 50 390.9 < 2.2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Intercept) D.PARK ## 4.3164849800 -0.0001058506
Y - Quantidade de veados infectados dentre os capturados
X - quantidade de vegetação aberta
plot(prop.veados ~ veados$openland,
ylim=c(0,1),
xlab = "Porcentagem de vegetação aberta",
ylab="Proporção de veados infectados")
nao.infectado <- veados$amostrado- veados$infectado y <- cbind(veados$infectado, nao.infectado) mod.veado <- glm(y ~ openland, data=veados, family="binomial") summary(mod.veado)
## ## Call: ## glm(formula = y ~ openland, family = "binomial", data = veados) ## ## Deviance Residuals: ## Min 1Q Median 3Q Max ## -3.8052 -1.1463 0.6707 1.3416 3.5483 ## ## Coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) 2.4820 0.1736 14.29 <2e-16 *** ## openland -6.0336 0.4432 -13.61 <2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## (Dispersion parameter for binomial family taken to be 1) ## ## Null deviance: 376.18 on 31 degrees of freedom ## Residual deviance: 121.30 on 30 degrees of freedom ## AIC: 204.65 ## ## Number of Fisher Scoring iterations: 5
anova(mod.veado, test = "Chisq")
## Analysis of Deviance Table ## ## Model: binomial, link: logit ## ## Response: y ## ## Terms added sequentially (first to last) ## ## ## Df Deviance Resid. Df Resid. Dev Pr(>Chi) ## NULL 31 376.18 ## openland 1 254.88 30 121.30 < 2.2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Intercept) openland ## 2.482035 -6.033583
Teste de Shapiro-Wilks
shapiro.test(arv3$alt)
## ## Shapiro-Wilk normality test ## ## data: arv3$alt ## W = 0.96262, p-value = 0.3609
shapiro.test(arv$alt[arv$sp == "sp1"])
## ## Shapiro-Wilk normality test ## ## data: arv$alt[arv$sp == "sp1"] ## W = 0.9906, p-value = 0.9591
shapiro.test(arv$alt[arv$sp == "sp2"])
## ## Shapiro-Wilk normality test ## ## data: arv$alt[arv$sp == "sp2"] ## W = 0.97453, p-value = 0.3502
Teste de Bartlett:
bartlett.test(arv3$alt~arv3$sp)
## ## Bartlett test of homogeneity of variances ## ## data: arv3$alt by arv3$sp ## Bartlett's K-squared = 1.1109, df = 2, p-value = 0.5738
Teste de Levene:
leveneTest(arv3$alt~arv3$sp)
## Levene's Test for Homogeneity of Variance (center = median) ## Df F value Pr(>F) ## group 2 0.512 0.605 ## 27
RESĂDUO = Y(observado) - Y(dados ajustados)
A premissa de normalidade na verdade recai sobre os resĂduos do modelo.
resid(anova.sps)
## 1 6 11 16 21 26 ## 0.8395548 -0.5666556 3.1201253 -1.1287012 5.8661704 -12.1406501 ## 31 36 41 46 51 56 ## 6.7500566 -4.5005793 0.1541789 1.6065001 -2.6382600 -3.1208224 ## 61 66 71 76 81 86 ## 3.3224765 -3.6806771 -0.3024985 6.3022481 -1.3052922 7.4459310 ## 91 96 101 106 111 116 ## -2.4542276 -3.5688778 9.2780228 1.1454976 5.9847796 -9.3347719 ## 121 126 131 136 141 146 ## 4.4754638 -2.2630830 -5.5808309 -6.9952785 0.7849457 2.5052549
#residuals(anova.sps) # o mesmo que a função de cima
plot do modelo# regressĂŁo: alt ~ nutri plot(alt.nutr)
# ANCOVA modelo aditivo: alt ~ nutri + sps plot(alt.sp.nutr)
# ANCOVA modelo interativo: alt ~ nutri * sps plot(alt.sp.nutr2)
# GLM poisson: atrop ~ dist.park plot(mod)
# GLM binoimial: prop.infectado ~ openland plot(mod.veado)